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Abstract 

For linear elastic problems, it is well-known that mesh generation dominates the 
total analysis time. Different types of methods have been proposed to directly 
or indirectly alleviate this burden associated with mesh generation. We review 
in this paper a subset of such methods centred on tighter coupling between 
computer aided design (CAD) and analysis (finite element or boundary element 
methods). We focus specifically on frameworks which rely on constructing a dis- 
cretisation directly from the functions used to describe the geometry of the object 
in CAD. Examples include B-spline subdivision surfaces, isogeometric analysis, 
NURBS-enhanced FEM and parametric-based implicit boundary definitions. We 
review recent advances in these methods and compare them to other paradigms 
which also aim at alleviating the burden of mesh generation in computational 
mechanics. 

Keywords: isogeometric analysis, review, mesh burden reduction, isogeometric 
boundary element methods, NURBS-enhanced FEM 



1. Introduction 

The finite element method is the most widely used numerical method in prac- 
tice and is underlined by fifty years of rigorous mathematical analysis. It offers 
results whose quality can be predicted a priori and controlled a posteriori for a 
range of problems of engineering relevance. The finite element method relies on 
the creation of a data structure known as "mesh" which is used to construct the 
approximation (of both the geometry of the domain and the field variables) and 
to perform numerical integration. The mesh is a set of elements of polyhedral 
shape (with straight or curved edges/faces) covering the domain. Traditionally, 
these elements have been simplex shapes (triangles and tetrahedra) or of quad- 
rangles/hexahedral shapes, but recently developed elements |14U] and numerical 
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Figure 1: Iterative design process. It is evident that such an iterative process is cumbersome 
and time consuming because a new, analysis-suitable, mesh must be generated for each new 
geometry. Several ideas have been advanced in the literature to try and decrease the mesh 
generation burden. 



integration techniques [103] allow the development of elements with an arbitrary 
number of edges/faces. 

For day to day industrial problems, initial analyses during early stage design 
are usually carried out using the assumption of linear elasticity. The geometry 
of the domain is most commonly generated using CAD software by a designer 
or engineer. A mesh must then be created to approximate this geometry and 
discretise the partial differential equation (PDE) governing the problem to allow 
analysis to be performed. The elements must fulfil a list of quality criteria 
related to their shape and grading throughout the domain, and this step is still 
difficult and far from being automated. Once the analysis results are available, 
it is typically necessary for the analyst to go back to the designer and propose 
changes to the design so that the part meet certain criteria related to, e.g. the 
maximum allowable stress or the maximum deflection in the component. In turn, 
a new mesh must be generated to approximate the new geometry and perform 
a new analysis. This iterative process must be performed as many times as 
necessary to converge to a suitable geometry. The idea is illustrated in Figure [T| 

It is evident that such an iterative process is cumbersome and time consum- 
ing because a new, analysis-suitable mesh must be generated for each new ge- 
ometry, which is not a simple task. Indeed, meshes of good quality elements 
are still today difficult to generate automatically, or are restricted to simplex 
(tetrahedral/triangular) elements, which are usually not robust (too stiff and 
not amenable to incompressible materials). Moreover, it may happen in indus- 
trial practice that the geometry is so complex that available mesh generators 



fail, or require significant human intervention. This has spurred a vast amount 
of research to enable analysts to generate quality meshes fast and robustly. An 
exciting avenue of investigation, followed by a large number of strong research 
groups worldwide consists in developing robust, fast, and automatic mesh gen- 
erators, in particular for hexahedral elements. Some of the most notable work 
performed in this area can be non-exhaustively summarised as 

Triangular/Tethedral-mesh generators Tetrahedral mesh generators are more 
flexible and robust for complex geometries compared to hexahedral mesh 
generators. The current meshing techniques can be roughly grouped into 
the following categories: advancing front |89l |88], octree-based methods 
|158j . Delaunay approaches [HH 1151] and mesh optimization techniques 
[5^ 03]. Some of the latest developments in this area include works on 
remeshing approach exploiting the natural anisotropy of most surfaces and 
the construction of high-quality isotropic tetrahedral meshes by using a 
variational principle [3i |4j. The reader is referred to the review paper of 
Alliez et al. [S| for details. The recent paper of Young et al. |159j is partic- 
ularly noteworthy as it is on the generation of meshes from medical images 
(image-based mesh generation). 

Quadrilateral/Hexahedral-mesh generators Hexahedral meshes are preferred 
from a numerical analysis perspective. However, it is more difficult to gen- 
erate, as explained in the survey of Shepherd and Johnson |135] . Hexa- 
hedral meshes can be generated directly, or converted from the tetrahe- 
dral meshes indirectly. Zhu [162] is one of the first ones using an ad- 
vancing front approach. Blacker and Stephenson [28] introduced a paving 
algorithm to provide a robust and efficient generator for all-quadrilateral 
meshes. Plasering is an extension of paving algorithm for hexahedral mesh 
in three dimensions. The reader is referred to the survey in 1998 by Owen 
[113] for unstructured meshes and refinement. Recent developments include 
the work of Zhang and Bajaj [160] who proposed an isosurface extraction 
method to automatically extract adaptive and quality hexahedral meshes 
from volume data, Ito et al. [75] who introduced an octree-based mesh 
generation method to create geometry-adapted unstructured hexahedral 
meshes automatically from triangulated surface models, Staten et al. [138] 
who developed an new algorithm to leverage the benefits of paving and 
plastering, and Sarrate and Huerta [122] who proposed an efficient algo- 
rithm for complex geometries based on a recursive domain decomposition. 

In parallel to these advances in meshing, several paradigms have been ad- 
vanced in the literature to try and decrease the reliance of numerical analysis on 
a mesh. The ideas behind these hinge upon the following goals summarised in 
Table [D 

Avoid altogether the generation of a mesh Meshfree methods (see Nguyen 
[104\ for a recent review and computer implementation details including 



an open source MATLAB code) relax the notion of neighbor relation- 
ships (connectivity) by simplifying the addition, suppression and move- 
ment of points/nodes. These methods came about in the seventies with 
the smoothed particle hydrodynamics (SPH) method proposed by Gingold 
and Monaghan [HS] and the element-free Galerkin method by Belytschko 
et al. Pl]. They are not as well understood mathematically as the finite 
element method, and most methods rely on the definition of a "mesh" , be 
it for numerical integration or generation of the point cloud. A rigorous 
quantification of the approximation property of meshfree methods, espe- 
cially with reference to the influence of the point distribution is an open 
problem, although some headway has been made, inter alia, in the area 
of point collocation methods developed by Davydov and Oanh |17j and 
maximum entropy interpolants by Rosolen et al. [120] . 

Element technology Meshes of tetrahedral elements are much more easily gen- 
erated than meshes of hexahedral elements (bricks) especially if user inter- 
vention is to be minimised. If robust and accurate tetrahedral elements 
could be formulated for most problems of engineering relevance, such fast 
mesh generators would become all the more useful. Several ideas have been 
proposed in the literature to avoid most of the limitations of (linear) tetra- 
hedral elements. The node-averaged tetrahedral elements of Bonet and 
Burton [29] , the locking free tetrahedral element of De Micheli and Mocellin 
[M] and the (cell, edge/face, node-based) smoothed finite element method 
(SFEM) of Liu et al. [SS] are examples of such attempts. Although these 
methods have other drawbacks associated with complications for multi- 
material problems and larger band widths, the edge-based smoothed finite 
element method was shown to yield superior results to that of the standard 
FEM for problems in linear elasticity [HS], visco-elastoplasticity |1U9] and 
geometrically nonlinear problems |1U8] . Polygonal/polyhedral elements are 
another alternative to ease the difficulties associated with generating qual- 
ity meshes made up of robust elements. They can serve as transition ele- 
ments during mesh generation, and therefore offer added flexibility. Work 
in this area includes the early work of Alwood and Cornes |6], Sukumar 
and Tabarraei |14Uj . Euler et al. [^I], Sukumar and Malsch |139] . Dai et 
al. [16], Mousavi et al. |101j and Natarajan et al. [103' . 

Use the geometry data provided by CAD directly in analysis From the 
description of the early-stage design process above, it is clear that the geo- 
metrical information provided by CAD is lost during the mesh generation 
process, because the smooth, arbitrary surfaces representing the boundary 
of the object are fitted as closely as possible during the mesh generation 
by piecewise polynomial surfaces (usually piecewise linear or quadratic) or 
the faces of the finite elements. The idea of subdivision surfaces B-spline 
finite elements of Cirak et al. [381 EZl EH] and isogeometric analysis of 
Hughes [To] is to use the geometrical description of the component to ap- 



proximate the field variables. More precisely, the (usually Non-Uniform 
Rational B-Spline - NURBS) description of the boundary of the domain 
is used to build shape functions for the approximation of the unknown 
fields on the domain. In this way, any modification of the geometry of 
the object is immediately translated into a modified approximation over 
the new geometry, without needing to reconstruct a polynomial fit of the 
boundary. The first set of such methods were developed in the context of 
the finite element method [70], where a domain mesh is required to con- 
struct the NURBS-shape functions, and, more recently, for the boundary 
element method (BEM) in the paper of Politis et al. |116] and Simpson 
et al. |136j . where mesh generation is completely avoided, bearing in mind 
the intrinsic limitations of the BEM. Advantages of the "isogeometric ap- 
proach" include the use of the exact geometry at all stages of the analysis 
and, more significantly, a significant reduction or circumvention of mesh 
generation. But recent research also illustrated several other advantages 
and difficulties of the IGA approach. It is the goal of the present paper 
to give an overview of some of these developments. Recently, the work of 
Moumnassi et al. |100j . based on the initial idea of Moes et al. [99] and 
Belytschko [25] showed that it was possible, using a marching method, to 
obtain an implicit domain definition based on multiple level set functions 
from arbitrary parametric surfaces provided from CAD data. This was a 
significant advance compared to previous work, as it enabled the treatment 
of corners and sharp edges (at the cost of a multiple-level set formulation). 
In some sense, this class of ambient space methods can be thought of as 
particular cases of the techniques described in the next paragraph. 

Use a structured mesh independent of the geometry An important dif- 
ficulty in the mesh generation process emanates from the requirement of 
the mesh to conform to the (usually complex) geometry of the domain. 
An alternative is to allow the geometry of the boundary of the domain 
to be non-conformally represented by the mesh. The boundary cuts the 
background mesh arbitrarily, and the latter is locally refined, usually using 
oct-tree based techniques - which also simplify parallelisation, in partic- 
ular using modern graphical processing units (GPUs). The method was 
recently revisited within the context of the extended finite element method 
[221 ESI [99]. 

Use methods only requiring a boundary discretisation This type of meth- 
ods take semi-analytical ways to decrease the dimension of the problem. 
The examples are Boundary Element Method initiated by Rizzo |119j and 
Scaled Boundary Finite Element Method by Song and Wolf |137j . Bound- 
ary element methods rely on fundamental solutions, which are not available 
for the non-linear problems and heterogeneous materials. On the other 
hand. Scale Boundary Finite Element Method can be extended to non- 
linear materials but not suitable for complex geometries where the scaling 



centre cannot be easily specified. 

In this paper, we focus on notions related to isogeometric analysis, in par- 
ticular to the isogeometric finite element method (IGAFEM), the isogeomet- 
ric boundary element method (IGABEM), the NURBS-enhanced finite element 
method (NEFEM) and some of the competing methods, in particular immersed/implicit 
boundary techniques. 

The paper is organised as follows: first, the technology that underpins CAGD 
and which forms the basis of IGA is outlined; next, the basic concepts of IGA 
and the fundamental differences to conventional implementations are described; 
recent developments in IGA are then described, emphasising the beneficial prop- 
erties over traditional approaches and finally, a comparative review with other 
methods which can also alleviate the mesh generation is given, followed by the 
comments on future directions and problems that IGA faces. 

2. Computer Aided Geometric Design Technology 

Isogeometric analysis (IGA) [70] has been successful at bringing two dis- 
parate research communities together: that of Computer Aided Geometric De- 
sign (CAGD) and numerical analysis. As argued above, converting CAD models 
into a form suitable for analysis is considered a laborious task, often dominating 
the design process. IGA offers a new direction, where, instead, a much more 
direct route is taken allowing analysis to proceed much more simply than before. 
This is achieved by using the data provided by CAD models directly rather than 
converting it through a preprocessing routine into a form suitable for analysis. 
IGA has found that the functions which are used to describe CAD surfaces are 
often directly amenable for analysis where, for example, properties such as Par- 
tition of Unity (PUM), linear independence, non- negativity and weak Kronecker 
delta property can be proven. 

To the authors' knowledge, the first work of integrating CAD and engineering 
analysis can be traced back to the paper by Kagan et al. [75], where B-splines 
were used as basis functions to represent both the geometry and the unknown 
fields. Following this idea, Cirak et al. [38J proposed a paradigm for thin-shell 
analysis, but used subdivision surfaces instead of B-splines. These ideas were 
formalised and generalised by the acclaimed work of Hughes et al. [70] based on 
NURBS, which has been extended to T-splines [12] recently. The recent work by 
Nguyen- Thanh |106l I105j and Wang |149] has incorporated PHT-splines. 

The current section is devoted to providing a non-exhaustive overview of the 
basic technology provided by the computer aided geometric design community, 
as this technology is central to isogeometric analysis. We attempt to perform 
this overview whilst keeping in mind that we wish to use the information pro- 
vided by CAD for analysis purposes. Table [2] gives a comparison of the technol- 
ogy reviewed here, namely B-splines, non-uniform rational B-splines (NURBS), 
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analysis-suitable T-splines, PHT-splines and subdivision surfaces. Other very 
recent techniques which are not reviewed here include RHT-splines. 



2.1. B-Splines 

B-spline geometry is a mapping from parametric space to physical space 
through a linear combination of B-spline basis functions, which are defined in 
parametric space, and the corresponding coefficients which are called control 
points because their physical meaning are a series of points scattered in physical 
space. B-splines in one dimension can be expressed as: 

n 

C(0 = 5^iV,,,(0B, (1) 

where C(^) denotes the physical curve we are interested, ^ the coordinate in 
parametric space, Bj the control points, Nj^p the B-spline basis functions of 
order p. 

B-splines possess the following properties 

The convex hull property The B-splines geometry is contained in the convex 
hull constructed by the control grid, which is a mesh interpolated by control 
points. 

The variation diminishing property No plane has more intersections with 
the curve than it has with the control grids. This property renders B-splines 
less oscillatory than Lagrangian polynomials. 

The transformation invariance property An affine transformation to a B- 
splines curve can be achieved by applying an affine transformation to the 
control points. 

The partition of unity property 

n 

1=1 

The non-negative property Each basis function is pointwise nonnegative over 
the domain. 

Ease of refinement and control of continuity It is possible to use the in- 
trinsic building blocks of NURBS for local refinement in the standard h 
and p way, but also by elevating the order of basis functions and continu- 



ity, simultaneously, also known as k refinement (see Section 3.2). 



Locally supported The basis functions are pointwise non-negative in p -|- 1 
knot spans. This property leads to a fact that the modification of a control 
point only infiuences the geometry locally. 
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Weak Kronecker delta property A weak Kronecker delta property means 
Nj(x.) = but Njijsij) 7^ 1, which is useful for enforcing boundary condi- 
tions in engineering analysis, because only the control points corresponding 
to boundaries need to be considered. 

2.2. Non- Uniform Rational B-Splines 

Non-uniform Rational B-Splines (NURBS) |115l 1124] are an extension of B- 
splines by introducing the weights for basis functions to represent a wider range 
of geometries such as conic sections. From a geometric point of view, NURBS 
in M*^ is a projective transformation of B-splines in M'^"'"^. NURBS inherit the 
properties of B-splines as mentioned above, but still have some drawbacks: 

Rational functions As NURBS are not polynomial functions, integrating them 
cannot be done exactly using Gaufi quadrature. 

Tensor product The parametric space and control points rely on a structured 
grid due to the tensor product property of NURBS, which increases the 
redundancy of the degrees of freedom. 

Continuity NURBS usually achieve only C° continuity between the patches. 

Geometry repair From a computational geometry point of view, NURBS based 
geometry always requires some level of repair due to gaps or overlaps of 
the various patches making up a geometry. 



The following subsections summarise a few of the recent advances to alleviate 
some of those shortcomings. 

2.3. Subdivision surfaces 

Subdivision surfaces, firstly proposed by Catmull and Clark [35], can possess 
a topology flexibility and ability to construct a watertight geometry. Subdivision 
surfaces use a recursive rule to interpolate or approximate a smooth surface. The 
methods can be divided into two classes: 

Interpolating Schemes When meshes are refined, new vertices are added with- 
out changing the vertices in the coarser mesh. Consequently, all of the 
vertices produced during subdivision will interpolate the surface. Interpo- 
lating schemes for quadrilateral meshes have been introduced by Kobbelt 
et al. [79], while Dyn [58] and Zorin et al. [163j described interpolating 
schemes for triangular meshes. In both cases the limit surfaces are but 
their curvatures do not exist. 



Approximation Schemes In approximation schemes, in contrast to interpo- 
lating schemes, the vertices in the coarse mesh need to be recomputed in 
refined mesh with adding the new vertices. The examples are the schemes 
proposed by Catmull and Clark [35] and Doo and Sabin |56], which are 
the first subdivision surface schemes and use quadrilateral meshes. Loop 
PU ] proposed a scheme based on triangular meshes. Theses approximation 
schemes can achieved continuous, except of some isolated points which 
are C^. After parametrisation, subdivision surfaces can be incorporated 
into IGA |38l [371139]. 

2.4- T-Splines 

T-splines were proposed by Sederberg et al. |126] to overcome the drawbacks of 
NURBS. T-spline control grids, also called T-meshes, permit T-junctions which 
are similar to the concept of the "hanging nodes" and oct/quad-tree meshes in 
FEM. In this way, lines of control points need not traverse the entire control grid. 
If a T-mesh is simply a rectangular grid with no T-junctions, the T-spline reduces 
to a B-spline. T-splines support many valuable operations within a consistent 
framework, such as local refinement, and the merging of several B-spline surfaces 
that have different knot vectors into a single gap-free model. T-splines form a 
subset of PB-splines, i.e. the basis functions are defined on local knot vectors. 
However, the local knot vector of T-splines can be inferred from a T-mesh, instead 
of being defined arbitrarily. A posteriori error estimation techniques for local h- 
refinement with T-splines have been very recently introduced by Dorfel et al. [FT] . 
Bezilevs et al. |12j have tested T-splines on some elementary two-dimensional 
and three-dimensional fiuid and structural analysis problems. 
The major drawback of T-splines is a significant implementation complexity. In 
addition, T-splines lose the property of linear independences, which is useful 
in engineering analysis. But an analysis-suitable T-splines P2], a subset of T- 
splines, can always satisfy this requirement. 

2.5. PHT-Splines 

PHT-splines (polynomial splines over hierarchical T-meshes) proposed by Deng 
et al. ^55] are constructed on the basis of T-splines and thus inherit their ben- 
eficial properties. PHT-splines use a hierarchical T-mesh which has a nested 
structure and different levels. PHT-splines define a cubic B-spline basis function 
in a T-mesh, which is level 0, and then modify the basis functions level by level. 
PHT-splines possess several main merits compared to T-splines: 

1. The basis functions of PHT-splines are polynomial. 

2. The local refinement is simpler compared to T-splines. 

3. The conversion between NURBS and PHT-splines is very fast. 

However, PHT-splines are only continuous, which are not sufficient to 
represent complex geometries. PHT-splines contain a closed set (linearly inde- 
pendent) of basis functions, which are important for combining them with FEM. 



3. Isogeometric analysis 



CAD models provide a complete description of the geometry of a problem, 
but for analysis to be carried out there is usually a need to carry out some level 
of preprocessing to arrive at a model suitable for analysis. This is also the case 
in isogeometric analysis (IGA), although the chief motive of this method is to 
simplify this preprocessing step. 

The difference between tradition practice and IGA however, is that the time 
and difficulty of this preprocessing step is significantly reduced. For example, 
in the case of IGA coupled with immersed and boundary element methods, the 
task of meshing is almost completely circumvented. This section aims to give 
an overview of the preprocessing steps required for IGA for various numerical 
methods and the differences and advantages over traditional techniques. 

By utilising the exact geometry and noting the advantageous properties of the 
basis functions that are used in CAD (see Section 2.2), several important devel- 
opments have been made in the computational mechanics community. Moreover, 
since the modifications required to convert a traditional numerical analysis code 
into one that utilises the same basis functions as CAD are often minimal, the 
method can be adopted quickly and applied to existing numerical methods. Sev- 
eral of these recent developments are highlighted in the present section where the 
benefits and drawbacks compared to traditional approaches are demonstrated. 



3.1. Motivation 

As discussed briefiy, the main purpose of isogeometric analysis (IGA), pre- 
sented by Hughes et al. [70] is to increase the ties between Computer Aided 
Design (CAD) data and the numerical method of choice for analysis. In the 
original paper, this method was the Galerkin Finite Element Method (FEM). 

For simplicity, let us take the example of linear elastic structural mechanics 
problems. In engineering practice, the geometry is described by CAD packages 
allowing the construction of very complex three-dimensional (3D) geometries, 
usually described by several NURBS (Non-Uniform Rational B-Splines) patches. 
This geometry then needs to be meshecj^ and analysed (e.g. by the FEM) to 
obtain an approximate distribution of the unknown fields within the structure. 
This analysis is used by the engineer to assess the design. If alterations to this 
design are required, they must then be communicated to the drafting team, and 
the geometry of the component suitably modified. A new mesh must then be 
regenerated by the analysis team to perform a second stress analysis, and the 
process is repeated as many times as is required for a suitable design to be 
obtained. 



^This meshing operation is inexact, and only as the mesh size goes to zero is the geometry 
exactly represented. This geometry approximation error however converges to zero faster than 
the approximation error on the unknown fields. 



Considering this design cycle, it becomes clear that tying the geometrical CAD 
data to the mesh in an automatic and direct way, where the mesh would be 
defined automatically from the CAD data is highly desirable. This desirable in- 
tegration is evident from the consolidation trends observed in the CAD / Analysis 
industry where CAD companies join forces with leaders in analysis, with the goal 
to streamline the design cycle and decrease lead time in computer aided design. 

In IGA the geometry can be represented exactly for analysis by direct use 
of CAD data through the use of NURBS-based approximations. This approach 
is being pursued widely both in the engineering community by improving the 
basis functions and seeking adaptive methods [831 [121 EZ] and applying the ideas 
to fields outwith mechanics, such as electromagnetics [33J. The idea is also 
being taken up by the applied mathematics community through the derivation 
of error estimates by, e.g. [2T]. This engineering and mathematics research has 
allowed a rather rapid development of the method through application-oriented 
investigations and sound mathematical analysis, respectively. 

3.2. Methods of refinement 

One refinement method is known as knot insertion, which is analogous to h- 
refinement in FEM and relies on splitting an element. However, the insertion 
of the existing values can reduce the continuity of basis function without pro- 
ducing new elements. The second is order elevation, which has similarities with 
p-refinement in classical FEM. Finally, fc-refinement has no analogue in FEM, 
because it allows increasing the order of the basis functions and the number of 
elements, simultaneously. The advantage of this method is that fewer degrees of 
freedom are introduced compared to p-refinement. |12] details the h-p-k- refine- 
ment schemes, Dorfel et al. [57] explains an adaptive /i-refinement, and Beirao 
da Veiga et al. [20] introduce error estimates for h-p-k- refinement schemes. 

3.3. A priori error estimation 

The difficulties of a priori error estimation in isogeometric analysis arises from 
two aspects: 

1. NURBS are not interpolatory. 

2. The support of NURBS basis functions are in general found to be larger 
than polynomial approximations of an equivalent order. 

The first difficulty is solved by pulling back the approximated function to the 
regular parametric space and noting that NURBS in M*^ is a projective transfor- 
mation of B-splines in M'^"'"^. The second difficulty is addressed by introducing 
the so-called "bent" Sobolev spaces [ISj . A meaningful approximation expression 
is derived to include not only the approximated function but also the gradient 
of the mapping, which does not change with the refinement. 

Finally, a significant result is that: the IGA result using p-order NURBS 
has the same order of convergence as the Finite element method using p-order 
polynomials. 



3.4- Bezier extraction 



To further integrate isogeometric analysis with existing FEM codes, Bezier 
extraction is utihsed for NURBS [31] and T-sphnes |123j to allow numerical in- 
tegration of smooth functions to be performed on C° Bezier elements. The idea 
of Bezier extraction is that localised NURBS or T-spline basis functions can be 
represented by a linear combination of the Bernstein polynomials. In addition 
to localising the support of basis functions into an element, Bezier extraction 
provides an element data structure suitable for analysis. That is, similar to La- 
grangian polynomials in traditional FEM, Bernstein polynomials do no change 
element by element. Furthermore, no intermediate parametric space is employed, 
hence the physical geometry can be mapped directly from parent elements. How- 
ever, it should be noted that Bezier extraction increases the degrees of freedom 
of the system. 

3.5. Isogeometric boundary element method 

As mentioned previously, the main bottleneck of IGA is the transition from 
the boundary representation to a domain representation suitable for analysis. 
An isogeometric approach using the framework of the boundary element method 
- coined the isogeometric boundary element method (IGABEM) - was proposed 
to solve this problem. The idea relies on the fact that both CAD models and 
boundary element methods rely on quantities defined entirely on the boundary. 
In the IGABEM, NURBS (and later, T-splines) basis functions are used to dis- 
cretise the boundary integral equation (BIE) for the surface geometry, boundary 
displacement and traction fields. The first implementations of this concept were 
pubhshed by Simpson et al. |136] and Politis et al. [116] for two dimensional 
linear elastostatic problems and exterior potential flow problems, respectively. 

Recently, IGABEM was extended to include T-splines for 3D linear elasto- 
static analysis in |125] with analysis on complex geometries demonstrated. More- 
over, Scott et al. |125] take advantage of a regularised form of BIE to simplify the 
implementation of singular integration inherent in boundary element methods. 
IGABEM inherits the benefits of both IGA and the classical boundary element 
method representing an attractive technique for early-stage design. 

IGABEM also promises important advantages for problems with infinite do- 
mains, such as acoustic and electromagnetic problems. Advantages are also an- 
ticipated in areas such as moving boundary problems (including fracture and 
contact) and inverse problems. Similar to classical BEM, however, IGABEM is 
not naturally adapted to non-hnear problems such as elasto-plasticity because 
of the requirement to evaluate a domain integral and thus, the application of 
the method is mainly restricted to linear problems. Compared to IGAFEM, the 
mathematical foundations of IGABEM require further investigation, particular 
with reference to methods that utilise collocation. 



3.6. Analysis- aw are geometric technologies 

IGA circumvents some of the difficulties associated with meshing but still have 
to overcome severe problems associated with CAD model quality. A main chal- 
lenge of IGA is that most of the currently available methods, to the exception 
of the isogeometric boundary element method |136] described in Section 3.5 are 
based on a volumetric model, whilst CAD models only provide a surface model. 
Furthermore, commercial CAD systems do not provide boundary representations 
of models in any general way to help construct full volume representations. The 
work to generate a trivariate geomertical representation is necessary. Interested 
readers can refer to [lOl EH [771 SH [H EHJ ^SM [EH] . 



3. 7. Applications 

3.7.1. Plates and shells 

There are two kinds of theories for plate and shell problems. One is the 
Kirchhoff-Love (KL) theory, which is suitable for thin plates or thin shells. In 
Kirchhoff-Love theory, rotations, which are the first derivatives of the displace- 
ments, are not interpolated directly, but through the displacement approxima- 
tion, hence leading to the requirement of a continuous representation approx- 
imation of the displacements, to yield rotations. To construct a element is 
not easy for polynomial basis functions. The second is Mindlin-Reissner theory 
for thick plates and shell problems, where shear-locking has to be dealt with 
carefully. In addition, the number of degrees of freedom is increased compared 
to the KL theory because rotations need additional interpolation. 

A rotation free analysis can be performed by IGA, i.e. the high order conti- 
nuity between elements can be achieved by only interpolating the displacement 
without approximating its derivatives. Another apparent advantage of IGA for 
shells is that it can capture the exact geometry to enhance the analysis accuracy. 
Although IGA for thin plates (shells) require a rotation boundary condition and 
continuity at finite angles, these constraints can be enforced by Lagrange multi- 
pliers or other similar methods. The application of IGA to plates and shells can 
be found in EOl EZl US 1142] , to which the interested reader is referred to for 
details. 

3.7.2. Fracture and damage 

The traditional Finite Element Methods have difficulties dealing with crack 
propagation problems because of the mesh regeneration adapted to the moving 
boundaries and the singularity near the crack tip. A great success in this area 
has been achieved by extended finite element methods proposed by Belytschko 
and Black [22]. The property of partition of unity of NURBS basis functions 
can lead to a combination of IGA and the extended finite element method. This 
combination is especially potent for complex geometries and curved cracks. First 
results were presented in the papers by De Luycker et al. |53] and Ghorashi et 



al. [HI]. Moreover, IGA also possesses advantages in implicit gradient damage 
models which necessitates a higher order fields. This problem has been circum- 
vented by the use of meshless methods as in, e.g. [71 1114[ H8]. Isogeometric finite 
elements allow for the construction of higher-order continuous field on complex 
domains, which was first demonstrated by Verhoosel et al. |146] . Other advances 
of the application of IGA in this area were presented by Borden et al. [52] using 
a phase field model, de Borst et al. for failure P9] and Verhoosel et al. jl45] for 
cohesive fracture. 

3.7.3. CFD and CEM 

The importance of the geometrical model in the numerical solution of Euler 
equations of gas dynamics in the presence of curved walls was identified in a 
Finite Volume context by [13 H] . The same problem when solving Euler equa- 
tions was identified in a discontinuous Galerkin framework in |l10] and a detailed 
study presented in |143j concluded that accurate results can only be obtained 
taking into account the curvature of the boundary. The benefit of using an exact 
boundary representation in the numerical solution of compressible fiow problems 
was illustrated in |13Ul I133 l|. showing that the error in aerodynamic quantities 
of interest can be reduced significantly when the approximated boundary rep- 
resentation of standard isoparametric finite elements is replaced by the NURBS 
boundary representation. The effect of accurate boundary representations in the 
numerical solution of viscous problems has been also demonstrated by [93l 180] , 
but the importance of the geometrical model is not exclusive of fiuid mechanics. 
Maxwell's equations are also very sensitive to an accurate geometric description. 
In |156] the error induced by isoparametric approximations of curvilinear geome- 
tries is studied. By solving the 3D Maxwell's equations in a sphere, the authors 
show that an exact description of the geometry reduces the error in one order of 
magnitude compared to isoparametric elements. In |129l I132[ 1131] several elec- 
tromagnetic scattering applications are used in order to demonstrate the benefits 
of using the exact boundary representation of the domain instead of the classical 
isoparametric description. In every example shown, IGA is more accurate than 
the standard isoparametric approach because of the geometrical accuracy. The 
application of IGA in this area can refer to the work of [21 [131 ttH US] in fiuids and 
turbulence, [H EHl HIl [13 [Ml [IZ] in fiuid structure interaction, and [Ml ITil] 
in electromagnetics. 

3.1.4. Contact 

Contact problems are still a challenging topic in computational mechanics 
because of the nonlinear problems arising from the boundary and the high sensi- 
tivity to the geometry accuracy, especially when sliding is involved. IGA satisfies 
the geometric accuracy requirement, however, the schemes for enforcing the con- 
tact boundary condition must be chosen carefully due to the fact that the control 
points are not located on the boundary. The methods fall into two classes. The 
first is the non-mortar method, which exerts contact boundary condition at the 



quadrature point. This scheme is also called knot-to-surface method in IGA 
because each quadrature point is represented by a parametric coordinate. On 
the other hand, mortar methods are a segment-to-segment technique and the 
boundary condition is imposed weakly with Lagrangian multipliers. The first 
work dealing with contact treatment in IGA includes [ISl EUl EB ESI IS21 E2] • 

3. 7. 5. Waves and vibrations 

One of the advantages of the IGA for vibration problems is that the NURBS 
basis functions are non-negative, thus the mass matrix is also non-negative. In 
addition, the finite element method results in a so-called acoustical branch and 
an optical branch. This branching is due to the fact that there are two distinct 
types of difference equations for the finite element: those corresponding to the 
end-point nodes at element boundaries, and those corresponding to mid-point 
nodes on element interiors. The acoustical branch corresponds to modes in which 
the neighboring end- and mid-point nodes oscillate in phase with each other, 
and the optical branch modes are the modes in which they are out of phase. 
Alternatively, the quadratic NURBS difference equations are all identical, and 
no such branching takes place. The reader is referred to |13] and [71] for details. 

3. 7. 6. Structural shape and topology optimisation 

In a sense, IGA seems very suitable for the shape optimisation because of the 
tight interaction with CAD. Indeed, whether the optimisation is done manually, 
by trial and error, or automatically, through optimisation algorithms, the fact 
that any change in the shape automatically results in a change in the approxi- 
mation of the field variables, without any additional mesh generation, makes the 
method attractive. Furthermore, the design velocity is accurate and the result 
does not require any smoothing process. Three of the relevant very recent papers 
in this area are |148l l36l 1128] . Nagy et al. |102] study the variational formulation 
of stress constraints in isogeometric design. 

A drawback of this method, which is a direct byproduct of its advantages, is 
that the geometry approximation is intrinsically tied to the field approximation. 
In such a formulation, the flexibility of modifying the field approximation inde- 
pendently of the geometry approximation is lost. Although no detailed study 
of the interrelationships between geometry error and field approximation error 
during structural optimisation is known to the authors, this flexibility may be 
important to ensure convergence to a globally optimal solution since it is not clear 
from the outset that the geometrical basis is optimal from an analysis standpoint. 

For topology optimisation, the initial method is cell-based or element-based, 
where the geometry can be represented by level sets |ll2i| or splines. The tradi- 
tional FEM encounters several main challenges. Firstly, in the fixed grid method, 
the optimal result depends severely on the design space which is constructed by 
the fixed grid. Secondly, postprocessing is required to transfer the end result to a 
CAD-suitable model, unless NURBS or other CAD-suitable approximations are 



used. Recent work in this direction was presented by [M]. NURBS-based IGA 
can avoid the difficulties of cell-based methods by using splines to represent the 
geometry. However, a trimmed-NURBS technique is necessary to deal with the 
topology changes |128[ 1127] . 

3.1.1. Locking 

It is well known that depending on the approximation properties of the numer- 
ical method, volumetric locking arises when dealing with (nearly) incompressible 
materials. A simple and effective method to overcome this difficulty is to use 
standard GauB integration for the volumetric strain and a lower integration or- 
der for the shear strain. However, this method is not suitable for anisotropic and 
large deformation problems. 

A more general way is B method proposed by Hughes [69], which is a three 
field method where the strain tensor is split into a deviatoric and dilatational 
parts, and the latter part is approximated by its linear projection into lower 
order space. Elguedj et al. [SD] have investigated B in isogeometric analysis for 
linear incompressible elasticity and F, a extension of B, for nonlinear elasticity, 
and concluded that IGA is more robust than the traditional j9-methods. 

Concerning shear locking, NURBS finite elements are not necessarily locking- 
free but this can be alleviated by a discrete shear gap method, as described in 

3.1.8. Phase field and gradient theory 

Both phase-field methodology and gradient theory typically lead to high order 
partial-differential equations. The Cahn-Hilliard equation governing phase field 
problems is a fourth-order differential equation. Consequently, the weak form 
necessitates continuous functions. To represent the size effect which cannot 
be explained by classical elasticity theory, the gradient theory was ffist proposed 
by Mindlin |96l |97]. Following this work, various versions of the method have 
been proposed but they share the common features that the strain gradient is 
introduced to the constitutive equation, yielding a higher-order equation than 
the classical elasticity theory. Traditional methods experience difficulties in con- 
structing high order fields, which can be solved by IGA readily. The application 
of IGA in phase-field and gradient elasticity problems can be found in [67] and 
[62] respectively. 

3.8. Challenges 

We have seen in the preceding sections that isogeometric analysis offers dis- 
tinct advantages over the standard CAD-Analysis-CAD design cycle. The main 
challenge of isogeometric analysis is to generate a three dimensional discretisation 
for the interior of the domain when knowing only the boundary representation, 
given by the CAD model. This volumetric parameterisation has been the topic 
of a number of recent research papers, and some of them are produced under 



Marie Curie Initial Training Networks in Europe. The reader is referred to the 
following work [H |9l] . 

In addition to the "mesh" generation difficulty mentioned in the previous 
paragraph, there are a few more stumbling blocks, which require further investi- 
gatiorj^ 

Integration Efficient numerical integration of the underlying approximation in 
isogeometric analysis is difficult since most of those are non-polynomial. 
This difficulty is common to enriched FEM and most meshless methods. 
The most recent advance known to the authors is the work of Hughes et 
al. [72], where the number of optimal quadrature points is chosen to be 
independent of the order of basis functions, and roughly equal to half the 
number of degrees of freedom. 

Mesh distortion The behaviour of isogeometric approximations for distorted 
meshes has been studied very recently by Lipton [83j, where it is pointed 
out that current finite element mesh quality metrics are too restrictive 
for IGA, i.e. that IGA is less sensitive to mesh distortion. Echter and 
Bischoff found that higher order continuity between individual patches 
can significantly improve the accuracy. 

Boundary conditions The lack of strong Kronecker delta property of isogeo- 
metric approximations (similar to many meshfree methods) leads to com- 
plications in the imposition of essential boundary conditions. A recent 
paper on this topic is [18]. Note that isogeometric approximations possess 
a "weak Kronecker Delta property." This property endows IGA an advan- 
tage over meshfree methods in the enforcement of boundary conditions. 

Error estimation and adaptivity As for any numerical methods, a priori and 
a posteriori error estimates and adaptivity are of prime importance for reli- 
able isogeometric analysis. A method based on h— refinement was proposed 
very recently by Dorfel et al. [57j. However, p— and k— (see [70] for a pre- 
cise definition) refinements are also possible, which offer great flexibility in 
the way the approximation space can be tailored to the unknown solution. 
One particularly interesting point is the measure of the adequacy of CAD- 
based approximations for analysis, in particular in the case of topology and 
structural optimisation. 

4. Other methods to suppress the mesh burden 

4.1. Meshfree Methods 

It was quickly realised that the finite element method was not best suited to 
all problems, in particular for large deformation problems and fragmentation. 



^The reader is referred to Table [Tl for a summary. 



The first meshfree method developed was the smoothed particle hydrodynamics 
(SPH) by Gingold and Monaghan |65] in the context of astrophysics problems. 

Meshfree/meshless methods also have advantages in the context of CAD to 
analysis transition, since they allow placing points/nodes with a lot more freedom 
than in the standard finite element method. This is because, in contrast to 
finite element methods, they do not employ elements in the construction of the 
approximation. Instead, a set of nodes, or centres, associated with a domain of 
influence are sufficient. The connectivity between the nodes, determined by the 
overlapping of these domains of influence can be defined more flexibly than in 
the FEM, and, in particular, can evolve in time. 

The basic idea behind most meshfree methods is thus to lift the strict connec- 
tivity requirements posed by the FEM. However, the arbitrariness in the node 
placement is relative since the quality of the approximation is known to be de- 
pendent on the geometrical location of the nodes. A recent and very interesting 
paper on this matter is that by Davydov and Oanh |17]. Based on these nodes 
and the size of their domain of influence, an approximation is constructed and 
can be enriched in various ways, to better represent the solution. The reader 
is encouraged to refer to the recent review [lU4j for details and an open-source 
CodeH 

Some of the most salient points to mention about meshfree methods may be 
briefly and non-exhaustively summarised as follows: 

• Mesh generation is vastly simplified, although work on the approximation 
properties of given point constellations is still ongoing. However, the influ- 
ence of the node placement is generally not as severe as in the FEM. 

• Arbitrary continuity and completeness can be achieved which can be useful 
for higher order continua or Kir chhoff- Love plate/shell formulations. 

• Discontinuities can be treated by modifying the weight functions, or through 
enrichment, similarly to partition of unity enrichment. See the paper of 
Nguyen et al. |1U4] for details. 

• Nodes can be added easily to locally adapt the approximation, as in the 
work of, e.g. Rabczuk and Belytschko |118] . 

Most meshfree methods may be recast as a particular case of the weighted 
residual method, with particular choices of the trial and test spaces. Depending 
on these choices and whether local or global weak forms are employed, a variety 
of methods can be constructed, with their respective advantages and drawbacks 
in terms of stability, accuracy, computational cost, convergence and robustness. 



■^The code, along with enriched isogeometric finite element and boundary element codes 
may be downloaded here : http://cmechcinicsos.users.sourceforge.net/ 



This is explained in detail in the papers by Atluri's group, see, e.g. [8]. 

IGA and meshfree methods share quite a number of similarities, in particu- 
lar the advantages of providing high order continuity. Compared to IGA, some 
meshfree methods (the reader can refer, for more details, to the well-known ar- 
ticles by Belytschko et al. [23], Rabczuk and Belytschko |117] . and Xiao and 
Belytschko |153j for a discussion on stability and accuracy of meshfree methods 
based on Eulerian and Lagrangian kernels) can be considered to somewhat alle- 
viate the difficulties associated with mesh distortion. On the other hand, IGA 
retains the exact geometry, and can more easily be incorporated into existing 
FEM analysis software, as was shown by recent releases of open-source imple- 
mentations by Vuong et al. |147] . de falco et al. [SO] and Nguyen- Vinh et al. 
|110] A combination with IGA and meshfree methods is possible by taking 
advantage of point based sphnes (PB-splines), as mentioned in [12]. 

4.2. NURBS-Enhanced Finite Element Method (NEFEM) 

In isogeometric analysis the whole domain is represented as a NURBS entity, 
and the solution of the boundary value problem is approximated with the same 
NURBS basis used for the description of the geometry. There are two main 
differences between isogeometric analysis and NEFEM. First, NEFEM considers 
the exact NURBS description for the boundary of the computational domain, the 
usual information provided by a CAD software. Secondly, NEFEM approximates 
the solution with a standard piecewise polynomial interpolation. Moreover, ev- 
ery interior element (i.e. elements not having an edge or face in contact with 
the NURBS boundary) can be defined and treated as a standard FE. There- 
fore, in the vast majority of the domain, interpolation and numerical integration 
are standard, preserving the computational efficiency of the classical FEM. Spe- 
cific numerical strategies for the interpolation and the numerical integration are 
needed only for those elements affected by the NURBS boundary representation. 

NEFEM was first presented for 2D domains in |129] . showing the advantages 
in front of classical isoparametric FEs using both continuous and discontinu- 
ous Galerkin formulations for the numerical solution of some test problems. It 
is important to remark that all the ideas presented in |129] are valid not only 
when the boundary of the domain is parameterized by NURBS, but for any 
piecewise boundary parametrization. The discussion is centered on NURBS 
boundary parametrization because they are the most extended technology in 
CAD. In [129J NEFEM was shown to be a powerful tool for solving compressible 
fiow problems governed by the Euler equations of gas dynamics. Several high- 
order FE methodologies for the treatment of curved boundaries are discussed 
and compared in |132] . including isoparametric FEM, Cartesian FEM, p-FEM 



"^See |http : / / cmechanicsos . users . sourcef orge . net/ to download the 3D open source 
MATLAB code for enriched isogeometric analysis. 



(a) 



Figure 2: (a) Knot lines of the NURBS surfaces defining the NASA almond, (b) surface 
triangulation, and (c) surface triangulation and knot lines 

and NEFEM. Numerical examples show that NEFEM is not only more accurate 
than FE methods with an approximate boundary representation, but also out- 
performs p-FEM with an exact boundary representation, showing the advantages 
of combining Cartesian approximations with an exact boundary representation. 
In [I3T] the extension of NEFEM to 3D domains is presented. We brief NEFEM 
in the following paragraph. 

Consider an open bounded domain Q whose boundary dQ, or a portion of it, 
is parameterized by NURBS curves in 2D or surfaces in 3D. In 3D every NURBS 
is assumed to be parameterized by 

S : [0,1? — > S([0,1]') cone 

A regular partition of the domain Q = [J^Qf. in elements is assumed, such 
that Qi f] Qj = 0, for i j. As usual in FE mesh generation codes, it is assumed 
that every curved boundary face belongs to a unique NURBS. That is, one 
element face can not be defined by portions of two, or more, different NURBS. 
Note however that the piecewise definition of NURBS is independent of the mesh 
discretization. Thus, the NURBS parametrization can change its definition inside 
one face, that is, FE edges do not need to belong to knot lines. This is a major 
advantage with respect to the isogeometric analysis. Moreover, different faces of 
an element can be defined by different NURBS surfaces. Figure [2] shows the knot 
lines of the surfaces that define the NASA almond and the surface triangulation. 
It can be observed that the spatial discretisation is independent of the piecewise 
NURBS surface parametrization. 

NEFEM possess the advantage of accurate representation of geometry and 
alleviate the difficulty of generating the interior isogeometric element. However, 
the difficulty of constructing high order continuous field still exist because of 
the use of polynomials. In addition, the control points are not related to the 
nodal parameters, hence it is not a geometry "interaction" language and loses 
the advantages in shape optimisation. 

4-3. Implicit Boundary Methods 

The implicit boundary method is a class of methods, including immersed 
boundary methods [98], fictitious domain [66j, embedded boundary [73], virtual 



boundary |121j and Cartesian grid methods |157j . These methods share the 
common features, that is, the meshes do not need to conform to the geometry. 
In addition, its advantage over IGA is that it allows a certain flexibility in the 
choice of basis functions, which may be different for the field variables and the 
geometry of the domain. However, it necessitates a "mesh geometry interaction" 
step in order to find the elements which are "touched" by the implicit bound- 
ary. The differences between various versions of implicit methods depend on how 
to construct the domain using a background mesh, explicitly or implicitly, the 
choice of integration scheme and the imposition of essential boundary conditions. 
Normally, the following three methods are used: Lagrange multipliers, penalty 
method, Nitsche's method (related to the family of augmented Lagrangian meth- 
ods) or modifying the basis functions to satisfy the constraint directly. 

Moumnassi et al. |100j propose an algorithm to define the implicit bound- 
ary of a domain from arbitrary parametric representations. More specifically, an 
automatic conversion from parametric surfaces to zero level sets for general un- 
structured meshes is presented by the authors. The use of parametric information 
can control geometrical errors at the boundaries which affect the convergence rate 
and solution accuracy, without changing the underlying fixed mesh for analysis. 
This algorithm helps integration with computer aided design and manufacturing 
systems. 

4.4- Smoothed Finite Element Methods 

Recently, novel finite element methods were born from the coupling of stabi- 
lized conforming nodal integration with the standard finite element method [86j . 
An overarching theory has been developed in the recent paper [84J. The main 
premise of this theory is the wish to achieve reliable results using lower order ele- 
ments, i.e. simple meshes (triangles, tetrahedra). SFEM retains the accuracy and 
inherit the advantages of triangular and tetrahedral meshes to represent complex 
geometries and can benefit directly from any advance in automatic remeshing. 

Furthermore, smoothed FEMs are a lot less sensitive to locking (volumetric 
and shear) as well as mesh distortion (because Jacobians are not required since 
no isoparametric mapping is used. In this sense, SFEMs are a way to improve the 
quality of the results obtained by simplex elements, thereby significantly reducing 
the need for human-intervention in the generation of hexahedral meshes. 

A few salient features of the smoothed finite element method can be sum- 
marised as follows: 

• Reduced sensitivity to mesh distortion. 

• Softer response compared to linear simplex elements (triangular and tetra- 
hedral) thereby allowing the use of unstructured meshes, much easier to 
generate automatically than structured hexahedral meshes. 

• Possibility to obtain upper bounds of the solution both compared to the 
displacement-based FEM solution and to the exact solution. 



• A range of methods producing results comprised between quasi-equilibrium 
methods and displacement-based finite elements. 

• Possibility to obtain "ultra accurate solutions" with the "a" FEM (see e.g. 
HH and [Ml). 

The basic ideas of strain smoothing in FEM and XFEM have been reviewed 
recently in [30] and the reader is referred to this paper for more detail. 

4-5. Solution Structure Method 

Solution structure methods (SSM) |134j can be viewed as a subset of the im- 
plicit boundary methods, but the idea to enforce boundary conditions is different 
and is worth briefing in this section. The key idea of SSM is to use a Kantorovich 
method, i.e. the solution of a differential equation with homogeneous Dirichlet 
boundary conditions can be represented in the following form: 

u = uj^ (3) 

where u is the solution of the differential equation in which we are interested, u 
is a known function that is positive in the domain and vanishes on the boundary, 
and $ is some unknown function, which can be discretised by 

n 

^ = J2Nih (4) 
1=1 

where Nj denotes the basis functions, which can be polynomials, B-splines, radial 
basis functions. In SSM, u is taken as the Euclidean distance function and the 
basis functions are B-slines. SSM can be thought of as a method that transforms 
any basis function to the ones which satisfy the essential boundary condition. 

In SSM, the distance function is approximated and the introduction of a 
weighting function increases the computational cost. Compared to IGA, SSM 
cannot guarantee geometric exactness because the distance function is normally 
approximated by the shape functions. 

4-6. Scaled Boundary Finite Element Method 

The Scaled Boundary Finite Element Method (SBFEM) of Song and Wolf 
|137] is a semi-analytical numerical method which does not require fundamental 
solutions. In SBFEM, the Cartesian coordinate system is transformed into a 
curvilinear coordinate system. The shape function is only dependent on the 
location of the surface. In two dimensional cases, the displacement is discretised 
by 



(5) 



where ^ is the radial coordinate and rj is the tangential coordinate, N{ri) is 
the shape function which only depends on t], hence only the boundary of the 
domain needs to be discretised. Substituting the above equation into the weak 
form of the governing equations leads to an ordinary differential equation in ^ 
which can be solved analytically. SBFEM retains the accuracy and inherits the 
advantages of the FEM and is thus suitable for non-linear problems. It also shares 
some of its properties with the BEM, that is, only the boundary of the domain 
needs to be meshed. It offers great facility for infinite domains problems and 
crack propagation problems. One special advantage of the SBFEM for fracture 
problems is that the crack surface does not need to be meshed. However, the 
SBFEM is not available for complex geometries in which cases the scaling center 
used to define the curvilinear coordinate system is not easily defined. 

5. Conclusions 

This paper critically reviewed some of the recent advances in attempting to 
integrate CAD and Analysis, with a special emphasis on isogeometric analysis 
and NURBS-enhanccd finite element methods. 

From this standpoint, it was said that isogeometric analysis, which shares 
a number of features with meshfree methods, although it has the important 
added advantage of conserving exact geometries, appears like a strong method 
to achieve tighter CAD- Analysis integration. For this to happen, however, a few 
points need to be addressed: 

• The basis functions must be the same as that used to describe the geometry, 
thus losing the fiexibility to take advantage of other types basis functions, 
which has particular use in some cases, such as in the spectral element 
method or in structural optimisation. 

• A general and efficient algorithm needs to be developed to construct the 
parameterisation of the interior of the domain, based on surface informa- 
tion. 

• Mixed isogeometric methods should be studied in order to successfully deal 
with instabilities and construct multi-field isogeometric elements. 

• CAD uses in general multi-patch NURBS boundary descriptions, which 
often lead to holes, gaps, trimmed surfaces, etc. T-splines do not have this 
drawback, and single patch T-splines are able to provide a description for 
most complicated geometries seamlessly. However, the implementation of 
T-splincs remains a challenge and it would most likely require the CAD 
community to embrace T-splines for T-spline based isogeometric analysis 
to have significant impact on current industrial practice. 

• IGA divide the elements in the parametric space, so it is difficult to capture 
heterogenous materials because they are scattered in the physical space. 



It seems likely that some combination of the various methods which were 
examined in this paper will provide winning strategies towards CAD-Analysis 
integration. The use of CAD-functions with Cartesian/Oct-tree meshes, for ex- 
ample, appears like an interesting avenue for investigation. 

Our open source codes for 3D enriched isogeometric analysis, extended and 
smoothed extended finite element method and extended meshfree methods can 
be downloaded from our Source Forge Page at |http : / /cmechanicsos . users . 
[sourcef orge . net/, 
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